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ABSTRACT 

We investigate the observable effects of feedback from Active Galactic Nuclei (AGN) on 
non-thermal components of the intracluster medium (ICM). We have modelled feedback from 
AGN in cosmological simulations with the adaptive mesh refinement code ENZO, investigating 
three types of feedback that are sometimes called quasar, jet and radio mode. Using a small 
set of galaxy clusters simulated at high resolution, we model the injection and evolution of 
Cosmic Rays, as well as their effects on the thermal plasma. By comparing, both, the profiles 
of thermal gas to observed profiles from the ACCEPT sample, and the secondary 7-ray emis- 
sion to the available upper limits from FERMI, we discuss how the combined analysis of these 
two observables can constrain the energetics and mechanisms of feedback models in clusters. 
Those modes of AGN feedback that provide a good match to X-ray observations, yield a 7-ray 
luminosity resulting from secondary cosmic rays that is about 10 times below the available 
upper limits from FERMI. Moreover, we investigate the injection of turbulent motions into the 
ICM from AGN, and the detectability of these motions via the analysis of line broadening of 
the Fe XXIII line. In the near future, deeper observations/upper-limits of non-thermal emis- 
sions from galaxy clusters will yield stringent constraints on the energetics and modes of AGN 
feedback, even at early cosmic epochs. 

Key words: galaxy: clusters, general - methods: numerical - intergalactic medium - large- 
scale structure of Universe 



1 INTRODUCTION 

Radiative cooling of gas in galaxy clusters is so efficient that 
most of the hot gas phase in their core ought to be removed 
on a time scale smaller than the lifetime of the system, produc- 
ing an inward mo tion of the cooling gas, a "cooling flow" (e.g. 
Fabian et al.lll984l) . However, dramatic cooling flows are not ob- 
served in real clusters, and additional non- gravitational heating 
mechanisms are likely to keep gas on a higher adiaba t (e.g. lKaisen 
1 199 it IPonman et ai1ll999l : iLlovd-Davies et alJl200oT) . It is largely 
agreed that active galactic nuclei (AGN) are a viable source of 
energy available for the self-regulation of galaxy clusters. Obser- 
vations show that the energy associated with AGN in clusters is 
in m ost cases sufficient to balan ce radiative losses in the ICM 
(e.g. IMcNamara & NulsenI |2007| . and references therein). How- 
ever, it is less clear how the energy is released from the compact 
(-Ckpc) region surrounding the central super-massive BH, to the ~ 
10—100 kpc cooling radi us. Very similar problems are encountered 
in elliptical galaxies (e.g. lSarazin & Whitdll987l : ICiotti & OstrikeJ 
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1 19971 : iBrighenti & Mathews! l2000l : ICiotti & Ostriked |2012|) . For 
clusters, an important issue is whether most of the energy in- 
put from AGN (or galactic winds) to the ICM has occurred 
much before the assembly of the clusters via "pre-heating" (e.g. 
i Bialek et al.l l200ll: iTqzzi & Normar] l200ll: IBrighenti & Mathews! 
l200ll : iMcCarthv et afl 2004 ), or at a low redshift within already 
formed clusters (e . g. | Binnev & T abor 1995i: IChurazov et"al I I2OO1I: 
Bruggen & Kaisej] l2002h . The first mechanism r equires a lower 
energ y budget, with total energies ^ 10 62 erg ( IMcCarthv et al.l 
1200 8h . while the s econd possibility require s energies in excess of 
~ 10 63 - 10 64 erg dMathews & GucfcOl lbh . However, there is am- 
ple evidence for strong AGN out flows in cool core clusters (e.g. 
IMcNamara & Nulsenll2007l|2012l) . 

Theoretical work suggests that the rea l evolution of heating 
in clu s ters might be a combin a tion of both (Brighenti & Mathews 
20061: IMcCarthv et al.l 120081 ; IVazzj l201ll : iDubois et al.l I2OI2F 
Short et alj2012h . 



Quasar-induced outflows at high-redshift, possibly follow- 
ing mergers of gas rich galaxies, have been observed (e.g . 
iNesvadba eUd1l2006l : I5unn et al.ll201ol : lBhattacharva et al.l 2008). 
At lower redshift, the mechanical work done by X-ray cavities 
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on the surrounding ICM may represent anot her viable mecha- 
nism for heating and mixing th e ICM (e.g. iDavid et all l200ll ; 
iMcNamara & Nulsenl I2OO7L l2012h . Additional mechanisms can 
modify the energy requirements of AGN feedback, by providing 
complementary heating/mixing on various scales. The heating from 
Alfven waves in cosmic rays (CR) enriched bubbles, and heating 
from Coulomb losses of CRs and the sur rounding thermal ICM is 
an additional intere s ting topic of research dLoewenstei n et al.ll99l |: 
Siiacki et all 120081: Icdafrancesco et aD l2004 iGuo & Oh] 12003 : 



Mathews & Gudl201 1 ilFuiita & Ohirall201 lh . 



Thermal conduction can somewhat reduce the energy bud- 
get that c entral AGN have to provide i n order to stem cooling 
flows (e.g.lBertschinger & Meiksirtfl98 6; Bregm an~& Davidl 19881 : 
iNaravan & Medvedevll2001r) . Magneto-rotational instabilities and 
heat-flux driven instabilities in the weak and anisotropic cluster 
magnetic field have been proposed to reduce the cooling of gas (e.g. 
Ouataert 2008; McCo urt et al.ll2oTlh . Finally, also major mergers 
have been suggested as a viable mechanisms to reduce the cooling 
catastro phe, even if the real efficiency of thi s mechanism is contro- 
versial dPoole et alj2006l : lBurns et alj2008h . 

Cosmological simulations with radiative cooling, star forma- 
tion and galactic winds are unable to reproduce the observed pro- 
files of gas temperature, m etallicity and entropy in the ICM (e.g. 
iKravtsov & Borganil2012l for a recent review). The most powerful 
feedback from AGN has been studied with cosmological simula- 
tions in the recent past. Several groups successfully implemented a 

treatment of thermal AGN feedback in cosmological GADGET sim- 
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ulations (e.g. Dalla Vecchia et al. 2004; Siiacki & Springel 2006; 

ISiiacki et al]|2007l : iMcCarthv et al.l |2010| ; iFabian et al.ll2010h . ~m 
the aforementioned papers, the growth of BHs at the centres 
of galaxies is followed using sink-particles, and the energy re- 
lease from each BH follows from the Bondi-Hoyle accretion rate 
("quasar" mode). The energy emitted by the AGN is released 
by heating up SPH particles surrounding the surrounding ICM. 
Similar methods have been i mplem ente d in AMR grid simu la- 
tions by ICattaneo & Tevssieil d2007h and iTevssier et al.l d201 lh in 

the RAMSES code. The creation of "bubbles" inflated by AGN 

I d 1 

during thei r "radio mode" was simul ated by IDalla Vecchia et al.l 
(2004) and Isijacki & Springell {2006) in GADGET. In this case, 
the energy released by the black hole (BH) is deposited within 
pairs of bubbles in the form of thermal energy or cosmic ray en- 
ergy, and exerts mechanical work PdVon the surrounding ICM 
while buoyantly rising in the cluster atmosphere. Other models 
of mechanical feedback from AGN have also been implemented, 
by assigning a "wind" drift velocity to gas particles surround- 
ing the BH, with velocities i n the range ~ 10 3 - 10 4 km/s 
( IDalla Vecchia & Schavel 120081 : IFabian et alj|2010h . iDubois et al] 
( l2010Ll201l[) implemented a scheme to follow bipolar kinetic out- 
flows ("jet" mode) from simulated AGN in RAMSES, m onitoring 
the growth of BHs using the same setup o flTevssier et alj J201 lh . In 
the framework of galaxy formation studies, run-time models of ra- 
diative and kinetic feedback from AG N in galaxy simulations have 
been developed bvlKim et al] d201l[) in ENZO simula tions, and by 
lOppenheimer & Pavel J2006T) and lGabor et al] ]201ll) in GADGET 
simulations. 

To date, however, little attention has been paid to the amount 
of non-thermal energy deposited in the ICM by the various feed- 
back mechanisms. This is an important aspect, since non- thermal 
emission from galaxy clusters can offer a complementary way 
of testing and falsifying feedback models . Once accelerated, CR 
hadrons can accumulate in galaxy clusters dBerezinskv et alj 19971) 



and produce a non-ther mal component that cou ld be detected by 
7-ray observations (e.g. lAckermann et al.l |2010T) . Secondary par- 
ticles are continuously injected into the ICM via proton-proton 
coll isions, possibly leading to detectable synchrotro n radiation 
(e.g lBlasi & Colafrancescoll999l ; lDolag & Enfflirj2000l) . The com- 
bined analysis of radio observations and 7-upper limits, however, 
presently suggests that most of the observed large-scale radio emis- 
sion in clusters cannot be due to secondary electrons , based on the- 
oretical estimates of the required total ener gy in CRs terunetti et al] 
l2007l:lDonnert et al.l2010bf ; lBrunettil2009h . and to the required val- 
ues of the magnetic field in the ICM, in conflict with observations 
( Brun etti et afl|2009l; IBonafede et al]|20ld: IPonnert et al] l2010iJ : 
iBonafede et al J2OI ll ; |jeltema & ProfumcfeOl lh . 

Therefore, the detection or lack of non-thermal emission 
from cluster centres may inform us about the energy budget 
of non-thermal particles, and of the history and modality of 
previous heating episodes in the ICM. In addition, CR parti- 
cles may have an important dynamical effect on the ICM (e.g. 
iRuszkowski & Ohl20l il ; lBrunetti & Lazarianl2011bh . and they can 
also a ffect the evolution of X-ray cavities powered by AGN jets 
(e.g. iMathews & Brighenti 2007; ISiiacki et al] I2OO8I: IGuo & Ohl 
2008 ; iMathews & GucfeOl l j) . 

In this work we study the observable non-thermal features re- 
lated to AGN feedback models in a cosmological framework. A 
few single-object simulations have been used to investig ate the 
role of CR feedback in stop ping cooling flows jGuo & Onl2008l : 
iFuiita & Ohiral l201lL l2012h . However, the acceleration (and re- 
acceleration) of CRs at shocks triggered by AGN, as well as 
at merger and accretion shocks have been neglected in previous 
works. 

To our knowledge, the present study is the first in which such 
detailed CR physics (e.g. particle acceleration, reduced fhermaliza- 
tion at the sub-shock, pressure feedback of CRs, effective adiabatic 
index of the baryon gas) as well as variety of AGN feedback models 
(quasar, jet and bubble modes) have been applied to cosmological 
simulations. 



2 NUMERICAL METHODS 

We have produced cosmological cluster simulations with the adap- 
tive mesh refinement code ENZO. On the basis of the public 1.5 
version of ENZO we have implemented our methods to model the 
evolution and feed back of CR particles injected at shock waves 
dVazza et al]|2012h , as well as our (simplified) implementations of 
energy release from AGN. 

ENZO is a grid and adaptive mesh refinement (AMR) code 
using the Piecewise Parabolic Method (PPM) t o solve the equa- 
tions of hydrodynamics, originally written by Brya n et al. I Jl99l) 
and developed by the Laboratory for Comput ational Astrophysics 
at the University o f California in San Diego jNorman et al]|2007l : 
ICollins et all201oh [ T l 

The detailed description of our modules for CR-physics and 
AGN feedback are presented in Sec. 2. 1-2. 1. In all runs in this paper, 
we adopted radiative cooling for a fully ionized H-He plasma with 
a constant metallicity of Z = 0.3Zp), and a cool ing function with 
a cut-off at T = 10 4 K JSarazin & Whitell 19871) , as in the public 
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version of ENZO, while the re-ionization background due to the UV 
radiation from early stars and AGN is modelled by keeping a gas 
te mperature floor (~ 2 • 10 4 K) in the redshift range 4 ^ z ^ 7 (as 
in lVazza et alfeoiol) . 

We did not include star formation and feedback through winds 
or supernovae in these runs, thereby reducing the complexity and 
memory usage of the code. Hydrodynamical simulations suggest 
that while supernovae are important to reproduce the observed 
metal distribution of the ICM, they do not have a signifi cant impact 
on th e thermal history of the ICM on large scales (e.g. IShort et all 
l20l2» . 

For the simulations presented here, we assumed a "concor- 
dance" ACDM cosmology withfio = 1-0, = 0.0441, Qdm = 
0.2139, fU = 0.742, Hubble parameter h = 0.72 and a normal- 
ization for the primordial density power spectrum as = 0.8. 



Kang & Jones 2007 
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2.1 Cosmic ray-physics 

The basic methods to model the injection, advection and pres- 
sure fe edback of CRs in our ENZO runs are explained in de- 
tail in IVazza et al.l fe012h . We assume that CRs are injected 
at shocks with an acceleration efficiency, r\{M), that only de- 
pends on the Mach number, M, which is given by diffusive 
shock accelera t ion (e .g. iBelll Il978l ; iBlandford & Ostriker 1 19781; 



Figure 1. Acceleration efficiency of CRs as a function of Mach number for 
different pre-existing ratios of P cr / Pg in the pre-shock region. From black 
to red, the different lines show the acceleration efficiency for P cr /Pg = 
0.0, = 0.06, = 0.12, 0.18, = 0.24 and = 0.3. 



ENZO runs, and compar es very well with th e shock finding methods 



Drurv & Voelkl 1 198ll; [Ellison et al.l Il995l; iKang & Jones! Il990t ; 

Malkov & Q'C Drurvll200ll : iKane & jSl200l ICaSEoij . P a P er on this sub -> ect <Vazza et alj | 2012 | ) , we add the treatment of 



we develop elsewhere (Vazzaetal. 2009). Compared to our first 



New CR energy is injected in the system by multiplying the energy 
flux through each shocked cell by the time step and the cell surface 
(at each AMR-level): 



E cr = r/(M) 



PuV s 



Ati 



(l) 



where p u is the pre-shock density, v s is the shock velocity, Ati and 
Aii are the time step and the spatial resolution at the AMR level- 
1, respectively. To ensure energy conservation the thermal energy 
in the post-shock region is reduced proportionally at run-time. The 
dynamics of the mixture of gas and CRs within each cell follows the 
total effective pressure of the two, P e fl = P g + P C r = p[(7~ l)e g + 
(7cr — l)e cr ], where 7 = 5/3, 7 cr = 4/3, e g is the gas energy 
density and e cr is the CR-energy density. The dynamical feedback 
of CR pressure is treated in the Riemann solver by updating the gas 
matter fluxes in the in 1-D sweeps along the coordinate axes and 



using the effective gamma factor (7 e ff 



(TjVtTcrfig) 



), for the 



computatio n of the local sound speed in cells. 

As in IVazza et alj d2012h . we use the relativistic value of 
7cr = 4/3 everywhere, which corresponds to the flattest possible 
momentum spectrum of CRs, through 7 cr = q/3 (with f(p) oc p~ q 
and q = 4 for j CI = 4/3). Fixing the value of 7 cr within the sim- 
ulated volume is an unvoidable assumption of the two-fluid model 
adopted here. However, once the CR-energy density is specified, 
the CR-pressure depends only weakly on t he spectral shape of 
f(p) and on the cut -off of the spectrum (e.g. Ijones &Kanslll993l ; 
ljubelgas et ai]|2008l) . 

To model the injection of CRs at each shocked cell, we devel- 
oped a run-time shock finder based on pressure jumps, similar to 
iRvu et all J200l) . At each time step, we flag cells with a negative 
3-D divergence, V • v < 0, and with concor dant local g radients of 
temperature and entropy, VS ■ VT > (e.g. lRvu et al.]|2003l) . The 
local Mach number is computed by inverting standard Rankine- 
Hugoniot jump conditions for gas pressure. This method can run 
in either fixed grid resolution or adaptive mesh refinement mode in 



several additional physical processes that are listed below. 



2.1.1 Shock re-acceleration 

We model at run-time the shock re-acceleration of CRs, by shocks 
running over a medium already enriched of CR energy by previous 
injections. This can be particularly relevant in the case of shocks 
caused by AGN feedback, where at late redshift the ICM is al- 



ready enriched with CRs ( P C r/-P g 



10" 



10 _1 at z < 0.5). 



According to the results in Kan g & Jonesl j2007l) the presence of 
CRs in the pre-shock region mimics an increased injection effi- 
ciency of CR-energy in the post-shock. This dynamic effect can 
be treated by using a different analytical function for rj(M), de- 
pendent on the ratio P CI /P g in the pre-shock. In our case, we cal- 
culate T](M) via a linear interpolation between the extreme cases 
of E CI j E g = (in which ca se the efficiency rj(Al) is identical to 
the one oflVazza et alj | 2012l) and P C i/P s = 0.3 (both taken from 
lKang& Jonesl 120071) . Figure Q] s hows the acceler ation efficiency, 
based on the interpolation of Ka ng & Jonesl d2007l) . as a function 

of the pressure ratio between CRs and thermal gas in the pre- 
shock region. The post- shock thermalizatio n is then reduced at run- 
time accordingly, as in IVazza et alj d2012l) . Based on our tests in 
IVazza et all (20 12). the run-time modelling of shock re-acceleration 
of CRs does not have a dramatic impact on the final distribution of 
CR-energy, that is on average increased by ~ 10 percent inside the 
virial radius. This is because the pressure ratio between CRs and 
thermal gas in always tiny inside clusters (Sec. 3.2), and therefore 
the impact of re-accelerated CRs on the final budget of CR-energy 
of the ICM is small. 



2.1.2 Hadronic and Coulomb losses 

Cosmic rays can lose energy via binary interactions with thermal 
particles of the ICM. This channel of energy exchange between 
thermal and relativistic particles in the ICM is important for the 
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high gas density (p/(nm p ) > 10 _2 cm~ 3 ) of cool cores. Rela- 
tivistic protons transfer energy to the thermal gas via Coulomb col- 
lisions with the ionized gas. They can also interact hadronically 
with the ambient thermal gas and produce mainly ir + , 7r~ and ir°, 
provided their kinetic energy exceeds the threshold of 282 MeV 
for the reaction. The neutral pions decay after a mean lifetime of 
w 9 ■ 1CT 17 s into 7-rays. To estimate the total energy transfer 
rate between CRs and thermal gas in both mechanisms, we need 
to determine the CR energy spectrum. Since this information is not 
readily available in the two-fluid model, we must assume an ap- 
proximate steady-state spectrum for the CR energy distribution. We 
fixed a spectral index of a — 2.5 for the particle energ)0 and com- 
puted the total Coulomb and hadronic loss r ates as a function of the 
ICM density and of E CI for each cell, as in lGuo & Ob] J2008I) : 



Tcoll = — Cc 



cm -3 erg • cm -3 



-1 -3 

erg • s • cm 



(2) 



where n c ~ n is the electron number density, and £ = 
7.51 ■ 10 -16 is the coefficient for all collisional energy loss terms. 
In hadronic collisions, only ~ 1/6 of the inelastic energy goes int o 
secondary electrons JColafrancesco et al.ll2004l : iGuo & Obi 2008). 
The energy-dominating region of CR electrons (7 ~ 10 2 ) will heat 
the ICM through Coulomb int eractions, plasma oscillations and ex- 
citation of Alfven waves (e.g. lGuo & Ohll2008h . Therefore we can 
assume that these secondary electrons lose most of their energy 
through thermalization and thus heat the ICM. Similar to Eq[2] the 
heating rate of the ICM through Coulomb and hadronic collisions 
can be computed as: 



I he 



E cr 



~ 3 ergs • cm~ 3 



-1 -3 
erg ■ s ■ cm , 



(3) 



where £ c = 2.63 ■ 10" 16 dGuo & Otfcoolj) . In our simulations 
with radiative cooling and AGN feedback, the rate of energy loss 
due to these collisions is extremely small, typically ~ 10~ 3 — 10~ 4 
of E CI or e g during the time step at each AMR level. This allows us 
to use a simple first-order integration to compute the energy losses 
of CRs (and the correspond ing gas heating rat e) at run-time. 

In our previous work jVazza et alj|2012l) . we did not include 
hadronic and Coulomb losses. In general, modelling this process 
at run-time decreases the CR energy by a factor ~ 10 within clus- 
ter cores, while yielding identical results for the remaining cluster 
volume compared to runs that neglect losses. 



BH growth using sink-par ticles Jsiiacki et alJl2007l ; lTevssier et al.l 
l201l| - |Martizzi et al.|[2012T) . we tuned the threshold value to excite 
AGN feedback to n m m rss 10~ 2 cm~ 3 . The use of this thresh- 
old is motivated by the fact that in our fiducial setup (peak res- 
olution of 25kpc/h per cell) the mass enclosed in a cell with 
n = 10- 1 - 10~ 2 cm -3 is « 10 10 - 1O U M . This is the typical 
gas mass surrounding BHs of Mbh ~ 10 8 — 5 ■ 1O 9 M0, which are 
com monly hosted inside the masses of galaxy clusters and grou ps 
(e.g. lSiiacki et alj2007l ; lTevssier et alj201ll ; lMartizzi et alj2012lrl 
In the following, we will refer to the cells exceeding this den- 
sity threshold and powering energy feedback as "AGN-cells". 



We have implemented three modes of AGN feedback: 
a "quasar" mode (i.e. thermal output of energy from AGN, 
Sec l2.2.lt . a "jet" mode (i.e. kinetic energy output from bipolar 
jets around AGN, Sec l2.2.2t and a "radio" mode (i.e. creation of 
buoyant bubbles in pressure equilibrium with the ICM, Sec |2.2.3l l. 

Once the feedback mode is specified, the only parameters that 
must be set are: a) the initial redshift for the start of AGN feedback 
(zagn) and b) the energy release of each single AGN-event, £"agn- 
In the case of "quasar" and "jet" modes this directly measures the 
energy we provide for each burst of either thermal or kinetic en- 
ergy, while for the "bubble" feedback this represents the estimated 
total energy released to the ICM by the creation of bubbles with 
internal pressure Pbb and volu me Vbb, Eagn « 3PbbVbb/2 (e.g. 
IScannapieco & B rtiggen 20081). 

Even if jets and radio bubbles are associated with the s ame 
type of AGN feedback (e.g. lMcNamara & Nulser]|200ll2012l) , in 
our study they are considered as alternative scenarios. This allows 
us to distinguish the CR-effects of buoyancy and impulsive kinetic 
feedback in a clearer way. 

In a preliminary set of tests (see the Appendix) we explored 
various recipes for the implementation of feedback modes in one 
reference cluster, before varying the efficiencies in the whole set of 
clusters. 

Here, however, we will discuss only the "fiducial" subset 
of parameters for which each implementation of feedback modes 
showed the best performance. While the range of spatial resolu- 
tion achieved in our runs is probably not sufficient to study specific 
small-scale features associated with each AGN mode (e.g. the mor- 
phology of jets or rising bubbles), our resolution and physical setup 
are suitable for studying the large-scale features of CRs in the ICM. 



2.2 Models of feedback from AGN 

In lVazzaldioi 1) we implemented a simple model of AGN feedback 
in ENZO, via injection of thermal energy at the opposite sides of the 
cooling region of galaxy clusters. Here, we explore more complex 
recipes of energy feedback between the cold gas and the surround- 
ing ICM, allowing also for a direct input of CR energy from AGN. 

At each time step, the identification of a suitable location 
of the central super-massive BH is based on the simple mea- 
sure of local gas over-density. First, we flag cells hosting a gas 
density exceeding a given threshold, n ^ iibh, and then we 
select as active "AGN-cells" only the maxima within cubic re- 
gions of size ~ 1 Mpc/h. Based on more detailed modelling of 

2 This choice of a corresponds to a Mach number of M = 3 at the particle 
injection. On average, this represents the typical injection spectra of par- 
ticles accelerated at powerful merger shocks in cosmological simulations, 
which ar e dominant source s of thermalization and CR injection within clus- 
ters (e.g. lVazza et alj201 lbl and references therein). 



2.2. 1 Thermal feedback from AGN 

A significant fraction of the energy emitted from AGN can ther- 
mally couple to the surrounding gas. One can define efficiencies 
such that the energy added in time At is 

ASagn.s = e r e c M B nc 2 At, (4) 

where e r ~ 0.1 is the bolometric radia tive efficiency for a 
Schwarzschild BH dShakura & Sunvaevll 19731) . Mbh, and E{ is the 
coupling efficiency with the thermal gas, which is usually assumed 

3 The choice of relying only on the gas density as a proxy may trigger also 
feedback from cold filaments connecting galaxies, where AGN feedback 
of the type considered here is unlikely. However, our simulations have suf- 
ficient high resolution only in galaxy clusters. This seems to exclude any 
spurious release of feedback energy from cold filaments, as indeed we find. 
To fully avoid this possibility, one would have to resort to more complex 
models (e.g. sink particles). 
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Figure 2. Time sequence of gas temperature (in units of [K] in the colour coding) for a slice of 5 X 5 comoving Mpc/h and width 25 kpc/h for run H5 with 
"quasar" feedback. The arrows point to the locations of recent "quasar" events in the volume. The "coarse" resolution of the first snapshots is due to the fact 
that we turn on AMR based on velocity jumps at z = 2. 



to be in the range ef « 0.05 — 0.15 in order to fit the observed 
Mbh vs <t v relation (e.g. Di Matteo et al.2005; Booth & Schaye 
2009). This thermal coupling between the AGN and the surround- 
ing ICM is usually called "quasar" feedback, which may be a mech- 
anism for pre-heating of the ICM JL"api et al . 2005; Mc Carthy et al.l 
|2004 120081: Isiiacki let alj 120091 : lLapi et al.l |2010|) . In addition, 
quasar-induced outflows at high-redshift, possibly following merg- 
ers of gas rich galax i es, have been observed in many cases (e.g . 
iNesvadba et alj|200d: Ehattacharya et alj|2008t iDunn et al.ll2010l) . 
Similar to McCarthy et alj d2010l) and lTevssier et alj 1 201 lh . when 
we detect cells with n ^ 7W n . we release the thermal energy, 
Sagn, adding to the total and internal gas energy inside cells at 
the highest available AMR level. 

This implementation of quasar feedback is only an approx- 
imation of the true physical processes at play, i.e. the launching 
of strong winds due to the radiation pressure of photons from the 
accretion disc. This is unavoidable, given that our best resolution 
is orders of magnitude larger than the theoretical accretion disc re- 
gion, and also given the difficulty of modelling the radiative transfer 
of photons from the accretion region. 

Our choice of a minimum gas density n m i n = 10~ 2 cm -3 
selects the typical environment of massive BHs within clusters and 
groups. Furthermore, we assume an energy output similar to most 
theoretical models, without actually measuring the accretion power 
of BHs at run-time. This is different from simulations where the 
mass growth of BHs is modelled using sink particles, which enables 
an accurate reconstruction of the BH matter accretion rate at run- 
time. This is not possible in our case. Our approach is only a first 
step to include the effects of AGN feedback in our version of the 
code, and it allows us a fast and efficient study of AGN modes. 



2.2.2 Kinetic feedback from AGN 

The innermost ICM can be affected by the injection of ki- 
netic energy through bipolar jets originating from the AGN (e.g. 
iBinnev & Ta bor 1995; Ciotti & O striker 2013). This energy can be 
thermalized by impacting on t he ICM afte r a short distance from t he 
clus ter centre (10 - 100 kpc.|Popell2009l) . lGaspari et alj d201 lbl!3) 
and iDubois et al.l d2010l 1201 1 ) recently simulated the mechanical 
coupling between purely kinetic jets from the AGN and the sur- 
rounding ICM in the cooling region. 

In this model each of the two jets is initialized as a pure in- 
put of kinetic energy density Eagn = l/2pvf ct , with velocity, 
Djet, pointing radially outwards from the cluster centre. Even if in 
our version of ENZO the launching direction of the jets is set by a 
random selection of coordinate axes, in this work we keep the jet 
axis fixed. In this way, velocity effects related to the direction of jet 
launching can easily be detected (Sec |3.4t . Every time a jet is gen- 
erated, we modify the gas velocity, the total and internal energy at 
the highest available AMR level in a pair of cells on opposite sides 
of the gas density peak. The width and the initial extension of the jet 
are set by the maximum resolution, which is ~ 10 ti mes larger than 
the bes t resolution available in "single-object" runs dGaspari et al.l 
2011b). At our resolution, the equivalent opening angle of each 
jet is ~ 30 — 40 degrees with respect to the cluster centre. This 
is an unavoidable drawback of the fact that in such cosmological 
runs achieving a much larger resolution is computationally very 
expensive. However, for the dynamical feedback of jets the most 
important quantity is the injected kinetic energy, which is similar 
to simulations at higher resolution, once the different mass load 
of the jets and launching velocity are considered. Also the initial 
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Figure 3. Time sequence of mass-weighted average gas temperature for a comoving volume of ~ (5 Mpc/h) 3 in the formation region of cluster H5. We 
show the results of the three re-simulations employing AGN feedback, always with an energy per event of -Eagn = 10 59 erg. Horizontal arrows suggest the 
location of interesting episodes of AGN feedback within the volume. 



velocity of our jets (which depends on the typical density of AGN- 
cells, through the total kinetic energy released in the AGN-burst) is 
~ 500 — 1000 km/s, about one order of magnitude lower than the 
typical velo city of jets in single- objects simulations at much higher 
resolution (Gaspari e t alj|201 lbl l3). 



2.2.3 Bubble feedback from AGN 

The creation of buoyant bubbles in the ICM inflated by jets re- 
quires a spatial resolution (~ 0.1 — 1 kpc) beyond what can be 
achieved in our cosmological runs. For this reason, we created al- 
ready formed evacuated bubbles aroun d cells hosting an AGN, sim- 
ilar to Scannapieco & Briiggen ( 2008). The bubbles are created in 
(approximate) pressure equilibrium with the surrounding ICM, by 
decreasing the gas density inside the cells by p' — <5bt>P, and cor- 
respondingly increasing the thermal energy in order to conserve 



the original gas pressure. During run-time, this modification is per- 
formed over a single time-step of the simulation, at the highest 
available AMR level. This obviously leads to a loss of gas mass 
in these cluster runs. The simulations show that on average this 
loss amounts to ~ 5 — 10 percent of the gas mass by the end of the 
simulation. 

The initial under-density of the bubbles is set by energy con- 
servation inside the bubbles, after that -Eagn is added to the gas 
and CR energy within the same volume. For the range of Eagn 
and n m i n adopted in this work, the under-density inside the bub- 
bles is typically 0.1 — 0.01. The bubbles are generated as a pair of 
cubic blocks of 2 3 cells, at the distance of 2 cells from the cluster 
centre. This is a very crude approximation and is expected to lead 
to expedient numerically mixing of the thermal energy. 

In principle, we can also allow for the presence of CRs in- 
jected by the AGN, parametrized by cf> CI — E C t,agn /Eagn of the 
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Table 1. List of the physical models adopted in our runs. Column 1: iden- 
tification name. C2: cooling. C3: energy per event. C4: details on feedback 
mode (with approximate values of heating temperature, jet velocity and ini- 
tial bubble under-density). 



ID 


cooling 


Sagn [erg] 


feedback mode 


cooling 


yes 





no 


quasars 


yes 


10 59 


thermal, T AGN ~ 5 ■ 10 7 K 


quasars2 


yes 


10 60 


thermal, T AG n = 5 ■ 10 s K" 


jets 


yes 


10 59 


kinetic, Vj ct ~ 750km/s 


bubbles 


yes 


10 59 


buoyant, 5j,bl ~ 0.05 



Table 2. Main parameters of the clusters simulated in this work. Column 
1: identification name. C2: total (gas+DM) mass at z = 0. C3: virial ra- 
dius (R v ). C4: dynamical state at Z = 0. The identification n ame of each 
cluster is ch oosen to be consistent with other works of ours IVazza et alj 
l201Qt[20T3) . The cluster parameters are referred to the simple non-radiative 
simulation of each object. 



ID 


A/tot [1O 14 M ] 


i?v[Mpc] 


dynamical state 


El 


11.2 


2.67 


post-merger 


E5A 


8.2 


2.39 


merging 


H5 


2.4 


1.70 


post-merger 


H10 


1.2 


1.20 


relaxed 



total injected energy, similar to lSiiackie t al. ( 2008). In this case the 
energy of the bubbles is conserved also by modi fying the effectiv e 
adiabatic index of the mixture of gas+CRs (as in lVazza et alj2012h . 
However, in the main body of this paper we only refer to purely 
thermal bubbles. We present a few results with gas+CRs inflated 
bubbles in the Appendix (the results, however, are not significantly 
different). 



3 RESULTS 

The results discussed in the main part of this paper have been 
produced considering the following physical mechanisms: a) pure 
cooling; b) cooling and thermal feedback by AGN from z ^ 4, 
with a fixed energy release per event of Eagn = 10 59 erg (or also 
10 60 erg in the case of El and E5A) ; c) cooling and kinetic en- 
ergy feedback by AGN from z ^ 4, with a fixed energy release 
per event of -Eagn = 10 59 erg; d) cooling and energy feedback 
from buoyant bubbles injected by AGN from z ^ 4, with an ap- 
proximate energy per vent of Eagn ~ 10 59 erg. With the above 
implementations, we re-simulated four clusters in the masses range 
10 14 M Q < M < 1O 15 M . The list of the most important param- 
eters of all "fiducial" models investigated in the article is given in 
Tab. 1. The first model (pure cooling) obviously produces strong 
cooling flows in all systems. In this work, this model is therefore 
regarded only as a standard reference to assess the role of each 
feedback model. 

In order to achieve the largest possible dynamical range in- 
side the volume where each cluster forms, the four clusters were 
re-simulated at high spatial and DM mass resolution starting from 
parent simulations at lower resolution and adding nested initial con- 
dition s of increasing spatial and mass resolution (e.g. lAbel et all 
Two levels of nested initial conditions were placed in cubic 
regions centered on the cluster centers. The box at the first level 
had the size of w 95 Mpc/h (with m dm « 5.4 • 10 9 M Q /h and 
constant spatial resolution of Ai w 425 kpc/h). The second box 
had a size of « 47.5 Mpc/h (with m dm ~ 6.7 ■ 10 s M Q /h and 
constant spatial resolution of Ai « 212 kpc/h). For every cluster 



run, we identified cubic regions with the size of ~ 6R V (where R v 
is the virial radius of clusters at z — 0, calculated in lower resolu- 
tion runs), and allowed for 3 additional levels of mesh refinement, 
achieving a peak spatial resolution of A « 25 kpc/h (in the fol- 
lowing, we will refer to this sub-volume as to the "AMR region"). 
From z — 30 (the initial redshift of the simulation) to z = 2, mesh 
refinement is triggered by gas or DM over-density. From z = 2 
an additional refi nement criterion based on 1-D velocity jumps 
dVazza et aljb009t) is switched on. This second AMR criterion is 
designed to capture shocks and intense turbulent motions in the 
ICM out to the clusters outskirts. Shocks and turbulence are ana- 
lyzed in the next Sections (Sec |3. 213.41 ). In addition, our test with 
ENZO have shown that this composite AMR criterion better cap- 
ture physical mixing motions in the ICM, and reduces the amount 
of numerical mixing, with important consequence s in the amount 
of cold gas within the cluster volume jVazzalboi ll) . However, this 
method is computationally much more expensive with respect to 
the gas/DM overdensity criteria alone, and we could afford to run 
it only in the "fiducial" AGN models (Tab. 1). 

In the Appendix we show the effects of different implementa- 
tions and energies budget in AGN feedback in the smallest of these 
systems,using the standard AMR criterion which works on gas/DM 
overdensity. In the main paper we will study the observable ther- 
mal and non-thermal features of the "fiducial" implementation of 
each feedback mode in our version of ENZO. The parameters of the 
small cluster sample resimulated with the different recipes of AGN 
feedback is given in Tab. 2. 

The evolution of gas temperature between z = 2 and z = 
in cluster H5 (quasar-mode) is shown in the panels of Fig(2] 
Blast waves triggered by AGN feedback (highlighted by arrows) 
in the thermal mode drive powerful shocks through the intra clus- 
ter medium, adding to the pattern of merger and accretion shock 
waves. We observe that in general only 2 — 3 AGN-cells are located 
inside forming clusters at a high redshift, while at a lower redshift 
only one AGN-cell is usually found inside the virial volume. The 
shocks triggered by AGN feedback are efficient in removing the 
cold gas phase at z > 1, while producing significant amounts of 
CRs in the forming cluster. Blast waves also cause the expulsion of 
a sig nificant fraction of the hot gas from the cluster volume, in line 
with lMcCarthv et alj feOlCh . 

Figure [3] shows a comparison of projected gas temperature 
across the entire AMR region of the same cluster, at three different 
cosmic epochs. While the large-scale morphology of the cluster is 
similar in all runs towards the end of the simulation, each feedback 
mode locally perturbs the ICM in a different way, driving power- 
ful outbursts (or buoyant bubbles) whenever strong cooling causes 
n > n min . 

In the quasar-mode, the final large-scale morphology of the 
cluster is also affected by the overall action of powerful thermal 
bursts, and presents a large amount of ejected hot gas outside of the 
cluster volume. In the jet-mode, prominent jets are launched along 
the horizontal axis, driving powerful shocks (M ~ 5 — 10) into the 
~ 10 6 — 10 7 K cooling gas as well as gas motions along their axis. 
The bubble-mode, on the other hand, has barely detectable effects 
on the large scale, and is globally similar to a simple pure-cooling 
run (not-shown). 

Fig.|4]presents a zoomed image of the gas density maps for a 
slice of 2 x 2Mpc 2 and depth 25 kpc/h for the final configuration 
of cluster E5A in runs employing AGN feedback. Jets and quasars 
are not active any more at z — in this run, and no clear features 
related to their activity can be seen in the innermost cluster region 
at this redshift. A pair of bubbles, on the other hand, has recently 



8 F. Vazza, M. Briiggen, C.Gheller 



been injected and can be seen as a couple of under-dense blobs on 
opposite sides of the cluster centre (shown by green circles in the 
image). Other images of active jets and bubbles can be seen in Fig. 
Al of the Appendix. 

We use phase diagrams for the cells inside the virial volume 
in order to characterize relevant differences. In Fig[5]we show the 
example of phase diagrams for an early cosmic epoch before the 
cluster formation, z — 2 (top panels), and after the formation of 
the cluster, z — 0.2 (lower panels) for the four different feedback 
modes. The most evident difference between the runs is the ab- 
sence of the "cooling flow" phase on the lower right part of the di- 
agram, in all runs employing AGN feedback. In the bubble-mode, 
the feedback from buoyant bubbles only affects the high-density 
regions around the cluster core. At high redshift, the phase diagram 
is very similar to the pure cooling run, with the exception of the 
selective removal of gas with n m j n > 10 -2 cm -3 cells, which are 
heated to high temperatures (> 5 ■ 10 7 K) by the injection of bub- 
bles. On the other hand, in the quasar and jet mode the action of 
feedback drastically alters the phase diagrams with respect to the 
pure cooling run. This is particularly evident in the quasar mode at 
2 = 2, where a "cloud" of high temperature cells (T > 10 7 K and 
n ~ 10 cm -3 ) is left after the passage of an AGN-driven blast 
wave. At z = 0.2, the "cooling" phase has almost disappeared in 
jet and quasar feedback, due to the efficient and volume-filling na- 
ture of the heating events driven by such powerful mechanisms. In 
the bubble run the phase diagram is similar to the pure-cooling run, 
with the exception of the n ~ 10 -2 cm -3 regime, where bubbles 
are injected almost continuously and over-cooling is balanced by 
mixing. This leaves a significant portion of "cold" (< 10 6 K) and 
dense cells in the simulated volume, characterized by large values 

Of Pcr/Pg. 

3.1 X-ray properties 

We compare our simulations with the collection of 
CHAN DRA cluster observations of Cavagno lo et al.1 
(2009), publicly available via the ACCEPT catalogue 
(http://www.pa.msu.edu/astro/MC2/accept). The catalogue 
consists of 241 clusters mostly located in the redshift range 
0.05 ^ z ^ 0.5 and with average central temperatures in the 
range 0.5 — 10 keV. For our comparison we selected a sub-sample 
of ~ 60 objects from the ACCEPT catalogue, with redshifts and 
average central temperatures compatible with our cluster dataset: 
z < 0.2 and 0.2 keV sC T < 5 keV inside 0.2i? 2 oo. For further 
details on th e ACCEPT catalogue an d set of observations, we refer 
the reader to Cavagn olo et al.1 J2009h . 

In Fig(6j we show the radial profile of reduced gas entropy 
Ekif/TsL, where K = T S i,/n 2/3 and T SL is the-spectroscopic 
like temperature. iRasia et alj 120051) . gas number density, n, and 
spectroscopic-like temperature for the different feedback modes. 
The simulated profiles are compared to the profiles of observed 
clusters, divided by the colour-coding into "cool-core" systems if 
their central entropy is < 30 keV cm 2 (in blue), or "non-cool-core" 
systems otherwise (red). 

As expected, none of the pure cooling runs can reproduce the 
observed profiles. In these runs, steep gas density profiles with cen- 
tral values of ~ 0.1 — 1 cm -3 and reduced central entropy be- 

4 We use the "reduced" entropy, rather than the usual entropy, K = 
T/n 2 / 3 , since this removes the dependence on the host cluster mass, 
lBorg anietalj2008l) . 



low ~ 10 cm 2 are found, at odds with observations and in agree- 
ment with standard radiative runs (e.g. lLi&Brv^nll2012l. for a re- 
cent review). Only cluster El is marginally within the range al- 
lowed by observations for density and entropy. However, consider- 
ing that this system has just gone through a strong major merger 
at (z ~ 0.25), it cannot be classified as a classical cool-core clus- 
ter, as suggested by its gas density and temperature profile. On the 
other hand, all runs employing AGN feedback yield a much bet- 
ter comparison with observations. Jet and bubble feedback modes 
produce roughly similar profiles for all clusters: a rather flat en- 
tropy profile within cluster cores (compatible with the h igh entropy 
floor of non-cool-core clusters in lCavagnolo et alj2009h . a shallow 
density profile and an almost isothermal temperature profile inside 
the central ~ Mpc from the centre. Less clear results are found 
in the case of quasar-feedback. While in clusters El, H5 and H10 
a reasonable match with observations is found with the "fiducial" 
energy budget of Z?agn = 10 59 erg per event, this is insufficient 
to quench the cooling catastrophe in cluster E5A, where a rather 
standard cooling flow takes place by the end of the simulation. By 
increasing the available energy per event by one order of magni- 
tude, £agn = 10 60 erg, the cooling flow is stopped also in clus- 
ter E5A (lower row of Fig(6}. Interestingly, the same is true for a 
re-simulation of cluster El using the same higher energy budget. 
However, in the case of the other smaller clusters the higher energy 
budget is too large and the thermal structure of both is destroyed, 
leading to a gas-poor cluster (see Appendix). 

This may simply suggest that different cluster masses must be 
characterized by different typical powers per event, mirroring the 
fact that a larger power per event is needed to balance the grav- 
ity an d the pressure of the cooling gaseous atmosphere (e.g. IPopel 
2009). However, it is interesting to notice that the self-regulating 
nature of AGN feedback (even with this simple implementation) 
yields very similar final profiles in the most massive object of our 
sample, even if the power per event is 10 times larger. In Fig[7]we 
show the total energy released by AGN in the formation region of 
cluster El, in the redshift range ^ z ^ 1. By the end of the sim- 
ulation, all feedback modes used a total amount of energy in the 
range ~ 2 — 9 • 10 61 erg, corresponding to ~ 0.01 — 0.09 of the 
total thermal energy of the cluster at z = 0. 

Interestingly, by the end of the simulation the higher power 
AGN mode (-Eagn ~ 10 60 erg) used about ~ 30 percent less en- 
ergy than the lower power run (Eagn = 10 59 erg), which is ~ 1/3 
of the total energy used in the kinetic feedback mode at a lower 
power. The total energy budget used in the bubble mode is simi- 
lar to the quasar mode at the same power per event. We note that 
in this system the use of AGN-feedback is significantly reduced 
for z > 0.6, due to the onset of large-scale mergers and power- 
ful shock heating within the cluster volume. However, episodes of 
AGN-feedback are present even at later epochs. Overall, a smaller 
amount of energy from AGN feedback is necessary to balance 
the catastrophic cooling with purely thermal feedback, or bubbles, 
compa red to kine tic feedback. This is consistent with analytical re- 
sults of Popg J2009L Fig. 1 ), where it is shown that a higher injection 
rate of kinetic energy, with respect to thermal energy, is necessary 
to balance the cooling flow for cluster masses > 10 M@. This is 
because in a massive cluster the critical momentum injection rate 
needed to overcome the pull of gravity and the surrounding gas 
pressure in the cluster core is more difficult to reach than the criti- 
cal injection rate of thermal energy required to overcome the cool- 
ing flow via thermal feedback. In addition, the energy release from 
jets close to the cluster core is very anisotropic, and it becomes 
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Figure 4. Gas density (top panels) for a slice of 2 X 2 Mpc/h 2 and width 25 kpc/h across re-simulations of cluster E5A at z = 0. The color coding is 
[p/Pcr bl (where p cr b is the baryon critical density). The green circles show the projected location of two bubbles previously injected in the cluster centre. 




Figure 5. Phase diagrams of cells within the virial volume of the re-simulations of cluster El at z = 2 (top panels) and at z = 0.2 (lower panels). The colour 
coding shows the cell-wise pressure ratio between CR and gas within cells, the thick line shows the average trend within each phase diagram. 



more isotropic (due to the driving of shocks) only at a distance of 
~ 100 - 200 kpc from the core. 

Overall, the X-ray properties of our clusters seem consistent 
with observed low-temperature non-cool-core clusters since the 
typical internal drop in temperature of cool-core clusters is not ob- 
served. On the other hand, in the epoch in which the AGN feedback 



is triggered in our simulations, the thermal structure of the inner- 
most cluster regions is not compatible with the appearance of typ- 
ical cool-cores either, because the internal drop in temperature is 
typically much stronger than what is observed. A lack of resolution 
compared to observations could be responsible for that. 

However, we find no evidence of bimodality in the distribution 
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Figure 7. Cumulative distribution of the AGN-energy released inside the 
comoving (4Mpc/fi) 3 of cluster El, from z = 2 to z = for 4 investi- 
gated feedback models. 



of clusters ent ropy or gas density in our re-simulations , at odds with 
observations dCavagnolo e7aiT2 009: Prat t et aljEjTcl :! Ecker t et al.l 
1201 ll) . This might be connected to our coarse sampling of the pa- 
rameters space of AGN feedback, or may instead call for the imple- 
mentation of more complex physics in our simulations. However, 
the investigated sample is too small to address this important issue 
in detail, and we will leave this for future work. 

We illustrate the effects of AGN feedback on the radial dis- 
tribution of baryon fraction by referring to the relevant example of 
cluster H5 at z = 0, shown in Fig(8] As expected, the onset of the 
cooling flow causes too strong a concentration of baryo ns inside the 
cooling radius, c learly at odds with observation (e.g. lEttori et all 
120091 : Sunl|2012l and references therein). On the other hand, all 
feedback runs produce a final distribution of gas mass which is 
more in agreement with observations , showing increasing profile s 
towards the cluster outer regions (e.g. lEttori et alj2(X)9tlSu32012h . 
Approaching -R200 the enclosed gas mass fraction is 90 — 100 per- 
cent of the cosmological baryon fraction. However, the different 
implementations of feedback produce significantly different and 
time-dependent features in the profile. In particular, the mechanical 
action of jets changes the shape of the enclosed baryon fraction in 
this low-mass system, due to a more recent feedback episode. Sim- 
ilar trends are found in the other systems, provided that the time- 
dependent action of jets of this power has a decreasing impact on 
the shape of the baryon fraction profile moving to larger masses. 
In the case of bubble-feedback, the artificial removal of gas when 
bubbles are created causes a small deficit of baryons inside the ra- 
dius. Since the enclosed baryon fraction is very similar (within a 
few percent at all radii) to the one obtained with thermal feedback, 
where the gas mass is conserved instead, we conclude that the loss 
of baryon mass due to our procedure of generating bubbles is not 
the leading mechanism of the observed trend. 
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Figure 8. Profiles of enclosed baryon fraction (normalized to the assumed 
cosmic value, /cosm), for the re-simulations of cluster H5 at z=0. 



3.2 Shocks and cosmic rays 

Blast waves triggered by explosive AGN feedback (in the quasar 
or jet mode) inject CRs at M > 3 shock waves in the ICM. This 
can happen before and after the formation of the cluster, as already 
shown in Fig[2] These shocks add to the underlying pattern of cos- 
mological shock waves driven by matter accretion, and to the bud- 
get of shock-accelerated CRs that are continuously injected into the 
ICM dRvu et alj|2003l ; IPfrommer et alfeood : ISkillman et ai]|2008l: 
IVazza et alj|2012l) . 

Figure [9] gives an example of the different patterns of shock 
waves and CR-energy in the re-simulations of one simulated clus- 
ter (we chose here the epoch of z = 1, when AGN feedback is still 
very active). The pure cooling run and the run employing bubble 
feedback show a pattern of shocks similar to non-radiative runs 
dVazza et"al]|2012h , with strong accretion shocks enveloping the 
cluster and a few weaker substructure shocks inside the cluster. The 
CR-energy in these runs is concentrated within the cluster and the 
filaments of matter being accreted onto it. The re-simulations with 
jets and quasars present a much wider distribution of CR-energy, 
as a result of previous episodes of gas and CR expulsion from the 
proto-cluster. In this case, shocks launched by mergers and accre- 
tions add to those previously driven by AGN bursts. In general, 
along the directions of powerful outflows of gas and CRs the loca- 
tion of accretion shock from smooth material is shifted up to several 
~ Mpc from the cluster center. This happens because the release 
of non-gravitational energy heats the ICM along the outflows, and 
increases the sound speed there. Very similar trends are found also 
in the other clusters. This implies that, in general, the impact of 
feedback on the distribution and amount of CR energy can be sig- 
nificant not only inside clusters, but also in large-scale filaments 
associated with the cluster. 

The differential volume distribution of shock Mach num- 
ber (measured as explained in Sec l2.lt for the same volume of 
(5 Mpc) 3 of Fig[9]is shown in Fig(T0j for all feedback modes. At 
the epoch of this output (z = 1), we find a systematic deficit of 
shocks in runs with quasar feedback, simply because shocks have 
travelled out of the reference volume, and the ICM is overall much 
hotter with respect to pure-cooling or bubbles feedback. The run 
with jets displays a distribution of Mach numbers closer to the pure- 
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Figure 6. Profiles of reduced gas entropy, gas numbe r density and spectro scopic-like temperature for the simulated clusters (thick lines with different styles) 
and for observed clusters in the ACCEPT catalogue iCavagnolo et al. 2009, the red lines are for non-cool-core clusters and the blue are for cool-core clusters). 
Each row shows the profiles of a different feedback mode. The last row shows the profiles for the quasar mode, with -Eagn = 10 60 erg per event (only the 
two most massive clusters are shown). 
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cooling case, but it shows a similar deficit of strong shocks inside 
the reference volume. 

In Figure[TT]we show the profiles of the pressure ratio P CT /P S 
for all clusters of the sample. As a reference, we also plot the av- 
erage pressure ratio obt ained with the sam ple of 7 non-radiative 
cluster runs presented in dVazza et alj|2012l) . 

In all feedback modes we measure a ratio P CI /P S ~ 0.1 in- 
side R v , within the range of what we already reported for non- 
radiative runs em ploying only CR injection at cosmological shocks 
dVazza et alj|2012h . However, large variations are found by com- 
paring pure-cooling or feedback runs. 

Cooling acts to increase this ratio by removing gas from its 
hottest phase at all radii. In addition, pure-cooling is also found to 
dramatically increase the loss-rate of secondary particles in the cen- 
tre, up to ~ 10 42 erg/s (as shown in Fig. A2 of the Appendix). The 
heating term from CRs is, however, not sufficient to balance the ra- 
diative losses in the innermost ICM, which are usually ~ 10 — 10 2 
times larger. Approaching the cluster core we measure P CI /P S ~ 
0.5 — 2 in the smallest objects, and ~ 1 — 10 in the most massive 
ones. Unfortunately, our sample of clusters is too small to allow 
us to infer a trend with the cluster mass. One possible reason for 
this trend can be, however, numerical: if the maximum spatial res- 
olution is fixed (as in this case) the cooling region of largest struc- 
tures is better resolved, usually producing a slightly larger gas den- 
sity peak and a significantly higher cooling rate. This trend in the 
cooling r ate has been reported in a number of w orks investigating 
SPH (e.g. Valdarnini 2002; iTornatore et al .120031) or Eulerian simu- 
lations (e.g. lBrvan& Norm an 1998; L i & Brvanll2012h . This trend 
with radius of th e pressure ratio is also similar to what is found al so 
with SPH codes jPinzke & PfrommeJl2010al ; lAleksic et alj|2012l) . 

We find that feedback acts to reduce the pressure ratio P CI /P g 
at all radii, this effect being larger inside R200 ■ This suggests that, 
although AGN feedback on average triggers the injection of ad- 
ditional CR-energy in the cluster core compared to simple radia- 
tive runs, its effect is that of reducing the dynamical role of CRs, 
because the same process also increases the thermal gas pressure. 
Although AGN feedback has the effect of changing the innermost 
shape of the profiles of pressure ratio in clusters H5 and H10 
(producing a profile more similar to the non-radiative case, with 
Pcr/Pg < 0.1 inside R200), the profiles of the pressure ratio in 
clusters El and E5A still show a maximum within the cluster core, 
with P cr /P g ~ 0.2 - 0.3. 

This internal shape is opposite to the case of the pure-radiative 
runs, and seems to be the "imprint" of the early cooling stage of 
both systems. We verified that indeed this profile is already in place 
at earliest epochs of cluster formation (z ~ 1 — 2) for both systems. 
It seems that once that the "cooling" profile had enough time to 
form, its imprint can hardly be erased by AGN feedback, especially 
in high-mass systems where, for a given AGN-power, it is more dif- 
ficult to produce an efficient outflow of the CRs previously accu- 
mulated. Similar profiles of the CRs to gas pressure ratio have been 
reported for sing le-object simulations of AGN feedba c k in c ooling 
flow clusters dGuo & Ohll200l ; iFuiita & Ohir all201 fl, |2012|) . sug- 
gesting that this may be a stable feature of such physical models. 
In order to study in detail the mass-dependence of P cr /P g in cos- 
mology, a much larger cluster sample is required. 

In general, it is not straightforward to immediately relate the 
power released by feedback (as in Fig[7]l to the injected CR-energy 
within the cluster volume. Shocks which follow the most recent 
event of jet or quasar feedback are surely responsible for the in- 
jection of new CR-energy inside the cluster, as shown in Fig(9] 
This can be seen also in the close correlation between the power 



per event of AGN feedback and the resulting P cr /P g pressure ra- 
tio inside clusters in the tests shown in Appendix (Fig. A3-A4). 
However, the long-term evolution of the simulated ICM and of its 
balance between CRs and gas energy depends also on other factors, 
such as the density-dependent secondary losses (Sec. 12.1.21 and the 
further evolution of matter accretion from the outer cluster region. 
Indeed, cooling and feedback do not only affect the cluster core, 
but also the density distribution in the accreted substruc tures, caus 



ing the release of different p attern of shock waves (e.g. lKang et al 
2007 ; Pfrommer et al]2007l) and turbulent motions (e.g. Valdarnini 



1201 II) . For this reason, it seems not possible to derive a simple pre- 
scription to relate the relase of feedback energy to the budget of 
CR-energy within clusters at low redshift, once that that all com- 
peting mechanisms are fully taken into account. 

We can compare these values to the available constraints from 
observations (see also next Sec l3.3l for a close comparison with 
7-ray observations). The present upper limits for a large sample 
of nearby gala xy clusters observed wi t h FERMI imply P CI /P S < 
0.05-0.1 (e.g. lAckermann et al.l20l"ol ; |jeltema & Profumol20lTh . 
with a poor dependence on the assumed spectra of CR. Also the 
limits on the presence of diffuse Mpc-scale radio emission in clus- 
ters can be used to constrain secondary el ectrons and thus the 
energy density of the primary CR protons dBrunetti et al J [20071 : 
iBrown et al 1201 ll) . In this case, the limits depend also on the clus- 
ter magnetic field strength and are complementary to those ob- 
tained from 7-rays. In the re levant case of an ave rage magnetic 
field in clusters of a ~ HG teonafedeet"ai1l2010l) , radio obser- 
vations of clusters with no Mpc-scale radio emission suggest that 
Pcr/Pg *S 0.05 (0.15) in the case of a spectrum with an index 
a = 2.5 (a = 3), whil e the limits are less stringent for smaller 
avera ge magnetic fields dBrunetti et alj|2007l ; jjejtema & Profumo 
I2OIII) . These limits usually refer to innermost ~ Mpc regions of 
clusters, where both the number density of thermal protons and the 
magnetic field are larger. At present no tight constraints are avail- 
able for the clusters outskirts where the CR contribution might be 
larger. 

According to Fig.Qj] it seems that all runs with un-balanced 
cooling produce a ratio P cr /P g that is incompatible with radio/7- 
ray observations. In runs with AGN feedback, we find a consistency 
with the above upper-limits in the two smallest size clusters, and a 
significantly larger pressure ratio in the two largest systems (fac- 
tor ~ 2 — 3 in the center). However, when the 7-emission of our 
clusters is directly compared to the available data from FERMI, the 
results of AGN feedback are below the upper limits. This is be- 
cause the 7-emission is more sensitive to the density (oc n 2 ) than 
to the pressure ratio between CRs and gas (oc P cr /P g ), and there- 
fore the net effect of AGN in our simulations is that of producing a 
lower 7-flux compared to upper limits, although producing a large 
Pcr/Pg in the center. In the case of radio upper-limits, the values of 
Pcr/Pg in our runs can be better reconciled with radio data if the 
average spectrum of CRs is a ~ 3 (or steeper). As we will discuss 
in Sec l3.3l this is actually the case of most of the CR-energy found 
in our cluster runs. 



3.3 7-ray flux 

CR-protons colliding with therm al protons of t he ICM pro- 
duce secondary 7-ray emission dDenniso 5 fl98Ch . Since once 
accelerated CR hadro ns can accumulate in galaxy clusters 
(Berezinskv et a il l 19971) and produce a long-lived 7-flux signature 



dPfrommer & EnfilirJ 2004 iMiniatil 12003); lAndo & Nagail I20Q: 
iDonnert et al. 201 Obi : Pinzke & Pfrommeill2010bl) . observations 
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Figure 11. Radial profiles of pressure ratio between CR and gas for our cluster sample, in the case of pure-cooli ng (top left panel) o r of different AGN-feedback 
models. The shaded grey area gives the result (average and ±l<r-scatter) of the non-radiative cluster sample of Vazz aet alj f2pT3) . 



this energy range can provide an important test-bed for feedback 
models. Indeed, feedback episodes producing similar X-ray fea- 
tures at z = might still contribute a different enrichment of CRs 
across a whole cluster lifetime, and therefore a different 7-flux from 
our simulated objects. In this Section we investigate the level of 
secondary emission resulting from the acceleration of CRs in the 
different modes of AGN feedback. 

In the two-fluid model used here, we cannot follow the spec- 
tral evolution of the CR population. Hence in order to compute the 
7-ray flux from hadronic collisions we must assume a spectrum 
for the CRs at the time of observation. Shocks which developed 
in the already formed dense and hot ICM (z < 1), largely domi- 
nate the energy budget of shock acceleration in clusters, compared 
to earlier times (z > 1). It is therefore reasonable that the resulting 
population of CRs at later epochs has a spectral index that is associ- 
ated with these Mach numbers. Using SPH simulations modelling 
the injection and energy evo lution of CRs at cosmological shocks, 
IPinzke & Pfrommej 0010b) recently derived an average a = 2.3 
injection spectrum for their sample of simulated clusters (for pro- 
tons in the range ~ 1 — 10 GcV). In our case, based on previous 
work with ENZO, we know that shocks that are mainly responsible 



for t he injection of CRs for z < 1 are in the range 2.5 ^ M ^ 3.5 
(e.g. lVazza et al.l201 lbl and references therein), corresponding to a 
particle injection spectrum of a ~ 3.3. This follows from the basic 
relation 



a = 2 



M 2 + 1 
M 2 — 1 ' 



(5) 



that applies to the case of non-radiative shocks. 

In the following we will base our computation of the 7-flux on 
the limiting cases of a = 2.5 (M = 3) and a — 3.3 (M = 2). 

Within the virial volume of clusters, the ir° decay is ex- 
pected to dominate the cluster emission above > 0.1 GeV with 
respect to other radiative mechanism (e.g. Invers e Compton of 
secondary electrons, non-thermal bremsstrahlung, Miniatil [20031 : 
Ando & Nagai 2008). We also neglect turbulent re-accelerat ion 
of CRs (e.g. iBrunetti et~al]|2007l ; iBrunetti & LazariarJl201lH) . 



We notice that in principle also the signal from Dark Matte r annihilation 
can contribute to the same energy range [1 - 10] GeV (e.g. Pinzke 
l201l|). However, this has not been observed yet (e.g. lAckermann 



20121) . and the effect obviously strongly depends on the assumed cross- 
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Figure 9. Maps of average CR-energy for a slice of depth 100 kpc/h in the 
four re-simulations of cluster H5 at z = 1. Each panel has a side 5 Mpc/h 
comoving. The colour coding is logio(EcR) (in arbitrary units), the over- 
laid colours show the location of M > 2 shocks for the same region. 
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Figure 10. Differential distribution of shock Mach numbers as a function 
of the AGN feedback mode, for the same simulated volume as in Fig[9] 



Following the forma lism of IPfrommer & Enfilinl d2004l) (see also 
iDonnert eta l, 2010b), we compute the emission from an asymp- 
totic spectrum of particle, with a spectral index a 7 ~ a. The 
threshold proton energy is set to E m i n = 1 GeV. 

Our results for each implementation of AGN are shown in 
FigHU where we show the absolute luminosity of simulated clus- 
ters in the [1 — 10] GeV energy range, (typically the most sensitive 
in terms of average differe ntial energy flux) and th e upper limits 
from FERMI observations o flAckermann et adl J2010h . as a function 
of the total cluster mass within -R500 ~ -R v ir/2. The results in this 
energy range do not depend much on the assumed spectrum of par- 
ticles, and the difference between a — 2.5 and a = 3.3 is only 
~ 30 percent in 7-flux (shown by different symbols). 

All pure cooling runs produce an excess of 7-emission com- 
pared to FERMI upper-limits, up to a factor ~ 10. This is mainly 
due to the large density of thermal targets for the proton-proton 
collision within the cooling radius. 

Runs with AGN feedback, on the other hand, are below 
present-day FERMI upper limits (with the exception of the quasar 
run of E5A, with .Eagn = 10 59 erg per event). The bubble- 
feedback seems to provide the lowest amount of 7-emission across 
the sample, with 7-ray flux about ~ 10 times below the FERMI 
upper limits within each mass bin. Jet and quasar feedback also 
produce 7-emission below the upper limits, but b y a smaller factor. 

Based on current upper-limits from FERMI dAckermann et al.l 

l2010h and the rather small sample of simulated clusters, it is not 
possible to reject any of the feedback models based on the predicted 
7-ray luminosity. However, this approach can put robust constraints 
on feedback models and energetics, complementary to the analysis 
of X-ray profiles. Indeed, while the radial profiles of gas entropy, 
temperature and density studied in Sec |3.1| cannot alone lead to a 
clear assessment of which feedback mode is more suitable to ex- 
plain observations, a ~ 5 — 10 times lower limit on 7-emission 
could already rule o u t the j et/quasar modes. 

iFuiita & Oh ira (2012) recently simulated the effects of in- 
jection, streaming and heating of CRs from the central AGN in 
the Perseus clusters, producing estimates of non-thermal radia- 
tion from secondary particles at several wavelengths. The 7-ray 
luminosities produced by our AGN runs for cluster H5 (which 
has a mass close to Perseu s) are of the same ord er of magnitude 
as in the fiducial model in IFuiita & Ohiral J2012h . ~ 10 41 erg/s. 
This is ~ 2 orders of magnitude below the FERMI upper limit at 
[1 — 10] GeV. Our AGN runs are similar to other cosmological 
simulations that estimate the 7-ray power from secondary pa rticles 
in clusters dDonnert et alj|20 10b; Pinzke & Pfrommer 2010af) , our 
results. Contrary to lPinzke & Pfrommej d2010al) , we do not excise 
dense overcooled clumps from the simulated volume and the con- 
tribution of shocks related to AGN activity is added to the contribu- 
tion of CRs from cosmological shocks. While a part of the gas mass 
locked into these dense clump would eventually lead to star forma- 
tion and produce cluster galaxies, these structures are the same that 
are responsible for the production of the un-realistic X-ray proper- 
ties (Sec |3.1| > and therefore their excision would just "mask" a dra- 
matic problem of the simulation. In addition, without the inclusion 
of star-formation processes in our simulation, it is difficult to estab- 
lish which part of these cold clumps would turn into galaxies, and 



section of DM, a still open topic, and we do not model this contribution in 
any way. Indeed, our simulations can be helpful in assessing the statistical 
emission floor of 7-emission expected from structure formation processes 
and AGN feedback. 
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therefore removing their contribution in 7-emission would be arbi- 
trary. The 7-ray power from the CRs accelerated at cosmologica l 
shocks only is slightly smaller than in IPinzke & PfrommedEoiOah . 
The fact that the total 7-ray power from our clusters is lower than 
in SPH calculations stems from the typicall lower innermost gas 
density resulting in Eulerian simulations, and also in the the differ- 
ent way in which SPH and grid methods model the acc retion of gas 
matter within cluster cores, as already investigated in lVazza et al.l 
ll2011bh. For the numerical reasons for that, we refer the reader to 
lAgertz et alj d2007» : IWadslev et all l200^ : lMitchell et afl d2009h . 

3.4 Turbulent motions 

Feedback from AGN may add turbulent energy to the already tur- 
bulen t state of the ICM. Unlike in our previous work IVazza et al.l 
with this set of simulations we can study turbulent motions 
from, both, the accretion of matter and AGN feedback (in its dif- 
ferent modes). 

In Figure[T3] we show a map of turbul ence reconstructe d with 
the same multi-scale filtering technique o flVazza et al.U2012h . con- 
sidering one run with jet-feedback (top panels) and in one with 
quasar-feedback (bottom panels). From left to right, the maps of 
gas density, total velocity and turbulent for a thin slice (25 kpc/h) 
through the cluster centre are shown. This system has an ongo- 
ing merger at this redshift, plus large-scale accretion of gas/DM 
matter along at least three filaments. These processes drive mo- 
tions with bulk velocities of ~ 1500 — 2000 km/s on ~ Mpc 
scales, and turbulent velocities of ~ 600 — 800 km/s on scales 
of ~ 300 kpc. In the run with jets an additional horizontal ve- 
locity pattern of ~ 600 — 700 km/s is present (indicated by an 
arrow within the panels). This is the outcome of a recent jet burst, 
launched ~ 0.1 Gyr ago. Even if the jet output dominates the ve- 
locity and the turbulent velocity field in the surroundings of the 
cluster core (< 100 — 200 kpc), making it anisotropic, the overall 
velocity and turbulence within the cluster volume are dominated by 
the large merger, with a pattern similar to the quasar- mode. 

In general, we find that feedback affect the turbulent veloc- 
ity field in our simulations in two ways. First, by driving turbu- 
lence locally to the AGN (e.g. along the jets or behind the wake of 
bubbles), in addition to the overall cluster velocity field. This is a 
general findings in simulations (e.g.[S cannapieco & Bri iggenl2008l ; 
iHeinz et alfeOlOhlGaspari et alj2011atlDubois et alfeOllT) . 

Second, by affecting the compactness of accreted substruc- 
tures and cluster satellites, thus affecting their ram-pressure- 
stripping and their efficiency in stirring the ICM. Overall, this 
makes the trend of turbulence in the different feedback modes non- 
trivial to identify, and close to the AGN highly time-dependent. 

This is confirmed also by the radial profiles of velocity and 
turbulence for the three most massive clusters (Fig|14t, While the 
velocity profiles at z = are quite similar (except that for the clus- 
ter core) in El, increasingly larger variations between modes are 
found in E5A and in H5. In all cases, the pure cooling run pro- 
duces the highest total an d turbulent velocity field. This is typica l 
of pure-cooling runs (e.g. iDubois et al J 1201 it iLi & Brvanll2012h . 
and a result of the transonic cooling flow. In addition, overcooling 
leads to more compact substructures, that ar e more efficient d rivers 
of tu rbulent motions duri ng their orbits (e.g. Val darninill201 ll) . 

ISanders et alj ( 1201 ll) recently computed upper limits from the 
maximally allowed amount of turbulence for a collection of 28 
nearby (mostly cool-core) galaxy clusters, by fitting a thermal 
multi-temperature spectrum to observed XMM-Newt on spectra. Fol- 
lowing the same procedure as in lVazza et al. we computed 
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Figure 12. Simulated 7-ray luminosity resulting from proton-proton for the 
1 — 10 GeV energy range, as a function of Afsoo- The large symbols are 
for a = 2.5 spectra and the small symbols are for a = 3.3 spectra. In the 
bottom panel (quasar mode), we additionally show with red symbols the 
result from the re-simulations of cluster El and E5A employing Bagn = 
10 60 erg per event. In each panel we co mpared the simulated d ata to the 
available upper limits from FERMI (from Acke rmann et alj201fj|) . 
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Figure 13. Maps for a 4 X 3 Mpc/h 2 slice of depth 25 kpc/h for cluster El at z ~ 0.2, simulated with jets (top panels) or quasar (bottom panels). The 
quantities shown are gas density (first column, log\o[p/ p cr ,b])> total gas velocity (central column, [km/s]) and the turbulent gas velocity (last column, [km/s]. 
Note that the turbulent velocities have been multiplied by 5 in order to use the same colour bar of the total velocity field). The horizontal arrows in the jet-run 
show the location of a recent jet event. 



the turbulent velocity field on ~ 30 kpc scales (which roughly cor- 
responds to the pr ojected volume available to the observations of 
ISanders et al.l201 lh for our cluster runs. Since the turbulent energy 
spectrum on the smallest spatial scales may be affected by numer - 
ical effects on the smallest scale (e.g. IPorter &"W oodward 19941) . 
the energy of the turbulent motions on scales ^ 30 kpc was calcu- 
lated analytically from the measured total power spectrum on larger 
scales assuming Kolmogorov scaling |f] 

In Fig[T5] we show the comparison between the simulated 
points at z — and the XMM-Newton limits for the galaxy clus- 
ters of the sample. As in I Vazza et all d2011al) . the values for the 
cores of simulated clusters are generally below the observed upper 
limits. However, the pure-cooling runs show a velocity slightly in 
excess of some observed objects. AGN feedback generally reduces 
turbulent velocities within the core (as shown above) and pushes 
the observed temperature from the core to higher values. Basically 
all observable turbulent velocities on the ~ 30 kpc scale are below 



upper limits (with the exception of the re-simulation of E5A using 
bubbles, in which the turbulent velocity is of the order of the 2 low- 
est available upper limits within the same temperature range). This 
statistics is much more time-dependent than the 7-ray fluxes, given 
the very small time scale associated with a change of the average 
innermost temperature and small-scale gas velocity, in response to 
AGN feedback events. In these runs, we verified that variations of 
factors ~ 2 — 5 in the average temperature, and of factors ~ 2 — 3 
in small-scale turbulence can easily be found by comparing epochs 
separated by ~ 100 Myr. Radiative cooling and dissipation of tur- 
bulent motions on the one hand, and violent heating and driving of 
outflows from the AGN on the other can cause such negative and 
positive fluctuations of the measured temperatures and velocities 
inside ~ 30 kpc. Taking this into account, we suggest that by using 
this technique it may be less straightforward to unveil the presence 
of specific AGN modes in the real observation. Nevertheless, this 
approach appears still useful to limit the available energy budget of 
AGN feedback. 



6 In this case the Kolmogorov slope maximises the possible contribution 
of unresolved structures of the velocity field, since in general a slightly 
steep er power-spectrum is mea sured in these simulations at the smallest 
scale jVazza et alj2 Q09. 201 la). Numerical dissipation makes it difficult to 
constrain the slope at the smallest scale with precision, and therefore our 
choice here already provides the most robust possible test against observed 
upper-limits of turbulence. 



A complementary viable option in the near future is repre- 
sented by X-ray spectroscopy of Doppler-broadened iron lines in 
the ICM, such as the Fe XXIII line at ~ 6.7 keV. S everal theoreti- 
cal works already tested t he impact of mergers (e.g. Sunvaev et al.l 
j2003t iDolag et all 120051 : iBruggen et ail 120051 : IVazzaet alJkOlCt 
Zhuravleva et al.l |2012l) . or AGN-feedback (e.g. Brtigg erTet al.l 
2005l : lHeinz et alfcOld) in shaping the Fe XXIII line with numeri- 
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Figure 15. Scaling between the average temperature and the mean velocity 
dispersion (extrapolated for < 30 kpc) for all cluster runs, and for the 
observational sample by Sand ers et alj feOllI) . The additional dotted lines 
show the dependence of the ICM sound speed with the temperature. 

cal simulations. In this case we can provide an overall view of the 
impact of various AGN modes on the shape of the Fe XXIII line in 
our simulated clusters. 

In Fig[T6]we show the simulated broadening of the iron line 
from a region of 100 x 100 kpc 2 centred on clusters El and E5A, 
considering the three lines of sights. We also included the effect 
of thermal broadening since feedback can simultaneously change 
the temperature of emitting gas along the line of sight. To normal- 
ize the effect of very bright emission in the cooling run, we nor- 
malize all lines to their integrated luminosity within the field of 
view. The general effect of AGN feedback is to broaden the emis- 
sion line, compared to the pure-cooling run, producing additional 
structures of high velocity tails (even if with very low associated 
luminosity) with a FWHM ~ 20 - 30 eV. The different feedback 
modes present the largest difference along the z-direction (which 
is the direction of the jets). Run El presents the largest difference, 
since jets were active more recently in this system. However, we 
find that in general disentangling the role of gravity-induced mo- 
tions and of non-gravitational processes is a very difficult task. 
Another problem is that the high-velocity materia l must entrain 
enough high-emissivity material to be detectable (Briigg en et al.l 
l200a iHeinz et alj[201ol) . In order to efficiently deal with "width- 
driven" or "separation-driven" features of the broadened emission, 
and to disentangle the different s ources of them, so phisticated tech- 
niques are being developed (e.g. IShang & Obll2012l and references 
therein). 



4 DISCUSSION 

In this work we simulated the evolution of thermal and non-thermal 
properties of a small sample of galaxy clusters with a customized 
versio n of the AMR code ENZO flNorman et all2007l : ICollins et all 
2010). T he basic features o f our two-fluid model have been pre- 
sented in IVazza et alj fe012l) , while the implementation of shock- 



re-acceleration of CRs, Coulomb and secondary losses and differ- 
ent modes of AGN feedback have been introduced in this work 
(Sec |2.1|2.2b . To our knowledge, this is the first time that non- 
thermal effects connected to cosmological shock waves and AGN 
feedback are studied with cosmological simulations. We tested the 
performance of a subset of recipes for AGN feedback (see, how- 
ever, the Appendix for a larger survey of models) against avail- 
able observations of thermal profiles derived from CHANDRA X- 
ray observations (Sec |3.1t , upper-limits of secondary 7-emission 
from FERMI observations (Sec l3.3b and upper-limits of turbulent 
motions from XMM-Newton observations (Sec l3.4t . Our main re- 
sult is that non-thermal observables allow us to assess the reliability 
of each feedback mode against available upper-limits, and is thus 
complementary to X-ray observations. For a more detailed discus- 
sion of the possible consequences of this method we refer the reader 
to Sec[5] 

Here we list the main limitations of our current numerical 
method, and the way in which our results can be affected or biased. 

Firstly, our recipe to identify AGN-like cells in the simulated 
cluster volume (Sec |2.2t relies on the heuristic choice of a threshold 
in gas density (n m i n = 10 -2 cm -3 ). This choice is based on a 
comparison with the typical gas density of environments hosting 
supe r-massive BHs, reported in numerical work by other groups 
(e.g. lSiiacki et al.f2007l ; lTevssier et al.ll201ll) . The lack of the more 
consistent use of sink-particles to model at run-time the growth and 
accretion rate of BHs can be addressed with future developments in 
ENZO 2.0 JCollins et al]|20irjh . 

Secondly, we re-simulated a few clusters with several imple- 
mentations of AGN feedback in order to compute the final budget 
of CR energy and their related 7-flux. These simulations are fairly 
expensive, and an extended survey of a large sample of clusters re- 
simulated with all variations of AGN feedback is currently beyond 
our means. For this reason, an investigation of the presence of a bi- 
modality in the distribution of cool-core and non-cool-core clusters 
in our dataset was not possible. This is left to future work. 

The recipe for radiative cooling is fairly idealized, and relies 
on the public implementation of equilibrium cooling in a fully- 
ionized H-He plasma with constant metallicity. Even though a more 
self-consistent inclusion of metal ejection from supernovae and 
galactic winds and line cooling can exacerbate th e cooling catas- 
trophe in radiative simulations I Duboi s et alj|201 lh . the overall en- 
ergy budget required from AGN feedback in our simulations (and 
its influence on the CR energy budget) should not be significantly 
affected. However, the inclusion of metals in our description of the 
ICM can provide an additional way of studying turbulent motions 
in the ICM, since successful implementations of feedback must 
yield the bimodality of metallic ity profiles observed in c ool-core 
and non-cool-core clusters (e.g. Leccardi & Molendi 2008). 

The neglect o f star format ion may artificially incre ase (by a 
~ 20-30 per cent. lSmifh et afl e.g. l201ll ; ICen & Chisaril e.g. l201lL 
and references therein) the mass of the hot gas phase in our runs 
here. At the same time, the inclusion of feedback from star for- 
mation and supernovae mig ht somewhat reduce the energy budget 
required from AGN (e.g. iTornatore et alj 120031 ; IValdarninil |2003t 
iNagai et al]|2007l ; Ishort et al.ll2012h . We will investigate this topic 
in the near future, using simulations with the effect of star forma- 
tion included. 

The inclusion of magnetic fields in these simulations is not 
expected to change the overall thermal structure of the ICM 
dDolag et alj 1 19991 : IXu et alj 120091 ; ICollins et alj l2010l) . It can, 
however, affect the mixing of jets and bub bles with the ICM 
dO'Neill & Jonesll201(]| : lMendvgral et alj2012h . 
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Figure 14. Radial profiles of total velocity (upper lines) and turbulent velocity (lower lines) for different re-simulations of El, E5A and H5. 
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Figure 16. Simulated Doppler broadening of FeXXIII line along three projection of cluster E 1 (top) and E5 A (bottom) at z = 0. Both thermal and velocity 
broadening and considered. The emission is extracted from a projected area of 0.1 X 0.1 Mpc/h 2 . In each model the emission is normalized to the total inside 
the area. 



A final important caveat concerns the spatial resolution of our 
runs. The minimum cell size of 25 kpc/h (almost uniform within 
the AMR region) is sufficient to capture sho ck waves and turbu - 
lent features related to cluster mergers (as in IVazzaetal.ll2011ah . 
However, it is barely enough to capture the interplay between the 
AGN region and the cluster core. Modelling the accre tion rate on 
the ce ntral galaxy at the centre of the cooling radius JLi & BrvatJ 
|2012|) . the small-scale interaction between jets (e.g. iHeinz et al.1 



l200dlMorsonv et al]20irj|;lGaspari et al]2012h or infla ted bubbles 
(e.g. lScannapieco & Bruggen||2008h Siiacki et al. i 2008), requires a 



resolution of the order of < kpc. Moreover, to resolve the tur- 
bulence excited by cluster mergers, sloshing and AGN-jets in the 
same simulation, one would need to cover scales ranging from 
Rvir ~ 3 Mpc down to the presumed scale of physical dissipa- 
tion at ~ 0.1 kpc, in a rather uniform way, for a range of scales of 
10 3 - 10 4 . 
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5 CONCLUSIONS 

What can be learned from this exercise of comparing thermal and 
non-thermal observables of simulated and observed clusters? We 
argue that the approach outlined here can put robust constraints on 
the energetics, duty cycle, mechanism and epochs of feedback from 
AGN. Even if with present data it is not yet possible to reject spe- 
cific implementations of feedback, it can be a powerful approach in 
the near future. The amount of CR injection after each AGN burst 
is assumed to be given by th e shock acceleration e fficiency of dif - 
fusive shock acceleration in dKang & Jonesll2007l : ICapriolill2012h . 
The other mechanisms of dynamical (e.g. advection with the fluid) 
or energy (e.g. secondary losses) evolution of CRs depend the dy- 
namics of the simulated ICM. All additional mechanisms of CR 
injection that we neglected here (e.g. direct injection from super- 
novae, galactic winds, magnetic reconnection) cannot but increase 
the budget of CRs in the ICM, even if not substantially. Therefore, it 
is likely that our limits of 7-emission from CRs are slightly low. We 
also r emark that other mechanisms of CRs, such has CRs diffusion 
(e.g. lHanasz & Leschll200d ; llubelgas et~aT1l2008l) . shou ld be negli- 
gible for the > 25 kpc/h scales of interest here (e.g. iBlasi et al.l 
120071 and references therein^. 

Therefore, in the near future a set of 7-observation signifi- 
cantly below these numerical estimates for jets/ AGN feedback may 
imply one of the following possibilities: 

• the real power per event for each AGN burst is substantially 
lower than what we assumed here (~ 10 44 — 10 45 erg/s) with a 
different duty cycle. However, such lar ge powers have been inferred 
from observation s in several cases (e.g. Mc Namara & Nulsen 2007; 
IStott et al]|2012L and references therein). 

• the acceleration efficiency of CRs at M < 10 shocks (which 
are typical in AGN bursts and cluster mergers) is lower than what 
we assumed here. Indeed, the details of particle acceleration for 
2 ^ M ^ 10 are not yet robustly constrained by theory due 
to the difficulty of modelling the large range of spatial and tem- 
poral sc^les_^nrolyed_in the diffusive acceleration at such shocks 
(e.g. iKang & Rvll2010h . More recently, several groups employ- 
ing particle-in-cells methods investigated additional acceleration 
mechanisms for particles at shocks (e.g. shock drift acceleration), 
suggesti ng the possibilities of a diff erent dependence on Mach 
number ( Gargate & S pitkovskvl l201lT) . In addition, there is now 
growing evidence that the non-linear diffusive shock acceleration 
model has to be slightly revised to fully explain recent data from 
supernovae JCaprioiil2012l) . 

• the actual quenching of catastrophic cooling in real clusters 
can happen over time through a mixture of "violent" (jets, quasar) 
and more "quiet" (e.g. bubbles) episodes. A physically motivated 
mixture of these different violent and quiet phases may in princi- 
ple quench the cooling catastrophe without injecting too large CR- 
energy in the ICM. Theoretical work, indeed, suggests that this is 



7 In lEnfilinetallfeOllI) it is argued that if CR can stream along a quasi ra- 
dial magnetic field much faster than the Alfven speed, the CR-energy den- 
sity in radio-quiet clusters would be greatly reduced. Howev er, it seems 
that such an hypothesis can be excluded on theoretical grou nds jAchterberel 
119811 : Schlickeiser 19891 1 1994 Farmer & Goldreichll2 004). and also based 



on a large number of Faraday Rotation lEnBlin & Vogtll2003l:l Murgia et al 
200 4IClarkj|2004lGuidetti et al.ll2008l:|Bonafede et alJl20ialVacca et al 



2010) and polarization I Bonafede et al. 2011) data from nearby clusters, 



that seem to exclude a systematic difference in the topology of the ICM 
magnetic fields between clusters with and without large-scale emission. 



a viable poss ibility (e.g. ISiiacki et alj|2007l : iMcCarthv et ai]|2008l: 
IShort et alj201^) . 

Based on the results presented in this paper, we argue that 
in the future a careful modelling of non-thermal observables will 
be important to improve our understanding of AGN feedback. De- 
spite the energy that is theoretically available for AGN feedback, 
the AGN cannot be arbitrarily impulsive or continuous because in 
both cases this will affect observable non-thermal phenomena (such 
as turbulence and the injection of CR-energy). In particular, given 
the extremely long time for Coulomb and hadronic losses of CR- 
protons in the ICM, the investigation of (lack of) 7-emission or 
radio-emission from secondary particles is a powerful tool to un- 
ravel feedback during the earliest phases of structure formation. 
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APPENDIX A: TESTS WITH VARYING PARAMETERS 

In this Appendix we show results of a more extended study of the 
variation of thermal and CR properties of one cluster run (H10), 
with a more systematic exploration of the parameter space of our 
implementation of AGN feedback. To save computational time, all 
tests in this case have been performed using the standard mesh re- 
finement strategy, where the grid is refined in regions where the 
gas/DM over-density is larger than 3 times the overdensity of the 
parent cell (the differences in the thermal properties inside -R200 
are very small compared to the AMR criteria used in the main ar- 
ticle, while in general the properties of accretion shocks and turbu- 
lent motions tend to be significantly affect ed by resolution is sues 
and numeric al dissipation, as discussed in Vaz za et al.l d2009h and 
IVazzal l l201lh ). 

Figure [Al] shows the projected gas density and a temperature 
cut for four re-simulations of this cluster for the pure cooling case 
(C) and for 3 feedback models with similar feedback energy budget 
(A2, K2 and B2). Despite the fact that the outer accretion shocks are 
poorly described because of the AMR strategy employed here, the 
innermost region where the AGN feedback plays a role are resolved 
similarly to the runs studied in the main article. 

Table Al lists all re-simulations of H10 performed, with de- 
tails on the adopted feedback model. In the ensemble of profiles 
shown in Fig |A2IA4l we show all profiles of thermal and CRs re- 
lated properties for all re-simulations. 

We compared one standard non-radiative run wi th CR feed- 
back only from cosmological shock acceleration (as in Vaz za et all 
I2012L run NRO), one run where also the Coulomb and secondary 
losses are modelled (NR1) and the pure cooling case with losses 
for CRs (C). 

Some results can be derived from this preliminary compari- 
son of feedback modes, in additon to those already discussed in 
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Table Al. List of the physical models adopted in our runs. Column 1 : identi- 
fication name. C2: cooling. C3: details of feedback model. C4: mnemonics. 



ID 


cooling 


feedback mode 


mnemonic 


NRO 


no 


no 


no cooling, no. CR losses 


NR1 


no 


no 


no cooling. ,CR losses 


C 


yes 


no 


cooling 


Al 


yes 


thermal, T AGN = 10 7 K 


quasar low 


A2 


yes 


thermal, T AGN = 10 8 K 


quasar med. 


A3 


yes 


thermal, T AGN = 10 9 K 


quasar high. 


AO 


yes 


thermal, T AGN = 10 8 K z ^ 1 


early quasar 


K2 


yes 


kinetic, Vj ct = 10 3 km/s 


jets high 


K4 


yes 


kinetic, Vfet = 10 2 km/s 


jets low 


B2 


yes 


buoyant, <5 bM = 0.1 


bubbles 


B4 


yes 


buoyant, <5bbl = 0.05 


bubbles 


BC 


yes 


buoyant, <5bbi = 0.05, </>cr = 10 


bubbles 



the main article, a) The adoption of a more standard AMR crite- 
rion still gives internal profiles and trends with AGN feedback very 
similar to the more expensive resimulations with AMR based also 
on velocity jumps, b) The thermal properties of the cluster are un- 
changed by the adoption of CR losses at run time (Fig A2). The 
amount of CR-energy within the innermost 100 kpc and the pro- 
duction rate of secondary particles (Eq[]} are however reduced by a 
factor ~ 10 due to CR losses. However, the inclusion of CR-losses 
is expected to reduce the cooling rate in radiative simulations, ow- 
ing to the larger gas density attained there. However, CR-losses are 
not able to stem the cooling flow, unless other heating mechanisms 
by CRs (e.g. heating from Alfven waves excited by the stream ing 
of CRs in a ~ 1 - 10/xG magnetic field lFuiita & Ohir dlioTlh are 
at work, c) In the case of quasar feedback, an energy per event of 
10 60 erg (i.e. 10 times larger than the fiducial model for H10 dis- 
cussed in the main paper) in a cluster of this mass creates a baryon 
poor system, and a unrealistically large ratio of CR to gas energy 
(~ 40 percent within the virial radius), d) For the quasar mode 
(Fig. A3), an energy per event of 10 58 erg (i.e. 10 times smaller 
than the fiducial model for H10 discussed in the main paper) can 
produce an acceptable match to the observed X-ray profiles, similar 
to what we reported for the two largest cluster masses (Sec. 13. lb . If 
AGN feedback is switched off too early (z = 1), a cooling catas- 
trophe develops only for a small region, 50 kpc. e) In the jet 
mode (Fig.A4), too slow jets (^ 10 2 km/s) are unable to quench 
the cooling flow, f) At this rather poor resolution, the investigated 
differences in the initial density contrast of bubbles (Fig. A4), or in 
their initial energy composition (i.e. thermal energy or CR-energy 
dominated) produce extremely similar result on the thermal and CR 
properties of the cluster at z = 0. 
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Figure Al. Top panels: slice of temperature (for a slab of 25 kpc/h through the cluster centre) for cluster H10 at z = 0, re-simulated with 4 different physical 
models. Each image has a side fa 5 Mpc/h. Bottom panels: slice of gas density for the same re-simulations of cluster H10 (note that the size of the images 
has been reduced (fs 2.5 Mpc/h), to highlight the features in the cluster core). The arrows show the location of buoyant bubbles in the B2 mode. 
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Figure A2. Radial profiles of gas entropy, gas density, gas temperature, pressure ratio P C r/Pg, baryon fraction and secondary injection for the re-simulations 
of cluster H10 employin g radiative cooling only , and non-radiative physics with and without CR losses. We additionally show as grey line the profile of the 
clusters in the sample of iCava gnolo et alj2009l) . 
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Figure A4. Same as in Fig |A2l but for re-simulations of H10 with different implementations of jet or bubble feedback. 



